Source code for hysop.backend.host.python.operator.derivative
# Copyright (c) HySoP 2011-2024
#
# This file is part of HySoP software.
# See "https://particle_methods.gricad-pages.univ-grenoble-alpes.fr/hysop-doc/"
# for further info.
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
# http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
from hysop.operator.base.derivative import (
FiniteDifferencesSpaceDerivativeBase,
SpectralSpaceDerivativeBase,
)
from hysop.backend.host.host_operator import HostOperator
from hysop.tools.decorators import debug
from hysop.core.graph.graph import op_apply
from hysop.numerics.stencil.stencil_generator import (
StencilGenerator,
CenteredStencilGenerator,
MPQ,
)
[docs]
class PythonSpectralSpaceDerivative(SpectralSpaceDerivativeBase, HostOperator):
"""
Compute a derivative of a scalar field in a given direction
using spectral methods.
"""
[docs]
def setup(self, work):
super().setup(work=work)
dA = self.dA
if self.scale_by_field:
aview = dA.compute_slices
self.scale = dA.sbuffer[aview]
else:
self.scale = dA
@op_apply
def apply(self, **kwds):
self.Ft()
self.compute_derivative()
self.Bt()
self.scale_derivative()
self.dFout.exchange_ghosts()
[docs]
def compute_derivative(self):
from hysop.constants import BoxBoundaryCondition
for nd_dkd in self.nd_dkds:
self.Ft.full_output_buffer[...] *= nd_dkd
[docs]
def scale_derivative(self):
out = self.Bt.output_buffer
scale = self.scale
if self.scale_by_field:
out[...] *= scale
elif self.scale_by_parameter:
out[...] *= scale()
elif self.scale_by_value:
out[...] *= scale
[docs]
class PythonFiniteDifferencesSpaceDerivative(
FiniteDifferencesSpaceDerivativeBase, HostOperator
):
"""
Compute a derivative of a scalar field in a given direction
using explicit finite differences.
"""
@debug
def __new__(cls, **kwds):
return super().__new__(cls, **kwds)
@debug
def __init__(self, **kwds):
"""
Initialize a FiniteDifferencesSpaceDerivative operator on the python backend.
See hysop.operator.base.derivative.FiniteDifferencesSpaceDerivativeBase for
more information.
Parameters
----------
kwds: dict, optional
Base class arguments.
"""
super().__init__(**kwds)
assert self.direction is not None
self.d = self.Fin.dim - 1 - self.direction
[docs]
def handle_method(self, method):
super().handle_method(method)
csg = CenteredStencilGenerator()
csg.configure(dtype=MPQ, dim=1)
stencil = csg.generate_exact_stencil(
derivative=self.directional_derivative, order=self.space_discretization
)
self.stencil = stencil
[docs]
@debug
def get_field_requirements(self):
stencil = self.stencil
G = max(stencil.L, stencil.R)
d = self.d
# set min_ghosts for input field
requirements = super().get_field_requirements()
for is_input, (field, td, req) in requirements.iter_requirements():
if field is self.Fin:
ghosts = req.min_ghosts.copy()
ghosts[d] = max(G, ghosts[d])
req.min_ghosts = ghosts
return requirements
[docs]
@debug
def discretize(self):
super().discretize()
d = self.d
stencil = self.stencil
assert self.dFin.has_unique_attribute("space_step")
stencil.replace_symbols(
{stencil.dx: self.dFin.discrete_fields()[0].space_step[d]}
)
self.factor = float(stencil.factor)
assert not stencil.is_symbolic()
[docs]
def setup(self, work):
super().setup(work=work)
dFin, dFout, dA = self.dFin, self.dFout, self.dA
iview = dFin.compute_slices
oview = dFout.compute_slices
self._in = dFin.sbuffer
self.out = dFout.sbuffer[oview]
self.iview = iview
if self.scale_by_field:
aview = dA.compute_slices
self.scale = dA.buffers[aview]
else:
self.scale = dA
@op_apply
def apply(self, **kwds):
"""Compute derivative."""
super().apply(**kwds)
stencil = self.stencil
_in, out, scale = self._in, self.out, self.scale
iview, d = self.iview, self.d
if self.is_inplace:
dtmp = self.dtmp
stencil.apply(a=_in, out=dtmp, axis=d, iview=iview)
out[...] = dtmp
else:
stencil.apply(a=_in, out=out, axis=d, iview=iview)
if self.scale_by_field:
out[...] *= scale
elif self.scale_by_parameter:
out[...] *= scale()
elif self.scale_by_value:
out[...] *= scale
self.dFout.exchange_ghosts()